{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[-0.017612, 14.053064],\n",
       "        [-0.752157,  6.53862 ],\n",
       "        [-1.322371,  7.152853],\n",
       "        [ 0.423363, 11.054677],\n",
       "        [ 0.667394, 12.741452]]),\n",
       " array([[-1.395634,  4.662541],\n",
       "        [ 0.406704,  7.067335],\n",
       "        [-2.46015 ,  6.866805],\n",
       "        [ 0.850433,  6.920334],\n",
       "        [ 1.176813,  3.16702 ]]))"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "\n",
    "#加载数据\n",
    "def load_data():\n",
    "    with open('简单分类数据.txt') as fr:\n",
    "        lines = fr.readlines()\n",
    "\n",
    "    x = np.empty((len(lines), 2), dtype=float)\n",
    "    y = np.empty(len(lines), dtype=int)\n",
    "\n",
    "    for i in range(len(lines)):\n",
    "        line = lines[i].strip().split('\\t')\n",
    "        x[i] = line[:2]\n",
    "        y[i] = line[2]\n",
    "\n",
    "    #以y区分两类x\n",
    "    x0 = x[y == -1]\n",
    "    x1 = x[y == 1]\n",
    "\n",
    "    return x0, x1\n",
    "\n",
    "\n",
    "x0, x1 = load_data()\n",
    "\n",
    "x0[:5], x1[:5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.08204613, 10.65488423]), array([-0.01613109,  2.95908379]))"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#求向量均值\n",
    "mu0 = x0.mean(axis=0)\n",
    "mu1 = x1.mean(axis=0)\n",
    "\n",
    "mu0, mu1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[0.9579532 , 0.54249108],\n",
       "        [0.54249108, 3.70851954]]),\n",
       " array([[0.9579532 , 0.54249108],\n",
       "        [0.54249108, 3.70851954]]))"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#求协方差矩阵\n",
    "def get_cov(x):\n",
    "    c = x.shape[1]\n",
    "    cov = np.zeros([c, c])\n",
    "\n",
    "    for i in range(c):\n",
    "        xi = x[:, i]\n",
    "\n",
    "        for j in range(c):\n",
    "            xj = x[:, j]\n",
    "\n",
    "            #协方差 = (x1 - x1.mean()) * (x2 - x2.mean()) / (n - 1)\n",
    "            xi = xi - xi.mean()\n",
    "            xj = xj - xj.mean()\n",
    "\n",
    "            var = xi * xj\n",
    "\n",
    "            n = len(var)\n",
    "\n",
    "            var = var.sum()\n",
    "            var /= (n - 1)\n",
    "\n",
    "            cov[i, j] = var\n",
    "\n",
    "    return cov\n",
    "\n",
    "\n",
    "#等价\n",
    "get_cov(x0), np.cov(x0, rowvar=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([[0.9579532 , 0.54249108],\n",
       "        [0.54249108, 3.70851954]]),\n",
       " array([[1.75563826, 2.04944407],\n",
       "        [2.04944407, 9.42458334]]))"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#协方差矩阵\n",
    "sigma0 = np.cov(x0, rowvar=False)\n",
    "sigma1 = np.cov(x1, rowvar=False)\n",
    "\n",
    "sigma0, sigma1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([1., 1.])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#初始化w\n",
    "w = np.ones(2)\n",
    "w"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(10.736930361702129, 10.736930361702129)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#投影后的均值 = 均值的投影\n",
    "w.dot(x0.T).mean(), w.dot(mu0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 58,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(5.751454888467757, 5.751454888467757)"
      ]
     },
     "execution_count": 58,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#投影后的方差可以使用投影前的协方差矩阵计算出来\n",
    "p0 = w.dot(x0.T)\n",
    "\n",
    "w.dot(sigma0).dot(w), ((p0 - p0.mean())**2).sum() / (len(p0) - 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.888466806593992"
      ]
     },
     "execution_count": 59,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#这个就是loss函数\n",
    "def get_J():\n",
    "    #投影后,两类的中心点距离最远\n",
    "    fenzi = np.power(w.dot(mu0) - w.dot(mu1), 2)\n",
    "    #投影后,两类的方差最小\n",
    "    fenmu = w.dot(sigma0).dot(w) + w.dot(sigma1).dot(w)\n",
    "\n",
    "    return fenzi / fenmu\n",
    "\n",
    "\n",
    "get_J()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([-1.38135267,  1.37410665])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#暴力求梯度法\n",
    "def gradient():\n",
    "    global w\n",
    "    upsilon = 1e-2\n",
    "\n",
    "    l1 = get_J()\n",
    "    w[0] += upsilon\n",
    "    l2 = get_J()\n",
    "    w[0] -= upsilon\n",
    "    gradient_w0 = (l2 - l1) / upsilon\n",
    "\n",
    "    l1 = get_J()\n",
    "    w[1] += upsilon\n",
    "    l2 = get_J()\n",
    "    w[1] -= upsilon\n",
    "    gradient_w1 = (l2 - l1) / upsilon\n",
    "\n",
    "    return np.array([gradient_w0, gradient_w1])\n",
    "\n",
    "\n",
    "gradient()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4.768051340025158 [-0.38135267  2.37410665]\n",
      "5.4261552382036315 [-1.72587658  1.90296174]\n",
      "5.426155141029574 [-1.72756646  1.90479546]\n",
      "5.426155157822327 [-1.72923578  1.90664113]\n",
      "5.426155174541138 [-1.73090184  1.9084832 ]\n"
     ]
    }
   ],
   "source": [
    "#训练\n",
    "for i in range(50):\n",
    "    w += gradient() * 1\n",
    "\n",
    "    if i % 10 == 0:\n",
    "        print(get_J(), w)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([26.8738025 , 13.79270217, 15.95380985, 20.38252513, 23.1817392 ]),\n",
       " array([11.32389099, 12.7950117 , 17.37850534, 11.74550404,  4.01073614]))"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#投影\n",
    "p0 = x0.dot(w)\n",
    "p1 = x1.dot(w)\n",
    "\n",
    "p0[:5], p1[:5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD4CAYAAADvsV2wAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAcmUlEQVR4nO3df4wc93nf8feHNFn3JAU2KcbWD/LoNoJRWbAVk2BsNHBk2IlpVrGSIAEkEK5guyBoRanTNoDlEEhSFUTRGnVrx0EVplL8g1cZKWzXgk3ZYtIGstHI1lEQZTmKbEoVJZqCdZKiSA4LyBSf/rFz4HK5s79mduY7M58XsLjb2bmd7+3tPTP7PM98RxGBmZm135q6B2BmZtVwwDcz6wgHfDOzjnDANzPrCAd8M7OOeFXdAxjl4osvjq1bt9Y9DDOzxjhy5MizEbFp2GNJB/ytW7eyvLxc9zDMzBpD0vG8x5zSMTPrCAd8M7OOmDjgS7pD0jOSHu5b9geSfijpwey2K+dnd0p6VNIxSbeUMXAzM5vONEf4nwF2Dln+nyPi6ux2aPBBSWuBPwLeC1wJ3CDpylkGa2Zms5s44EfEvcDzM2xjB3AsIh6PiJeBLwDXzfA8ZmZWQBk5/JslPZSlfF475PHLgKf67p/Ilg0laY+kZUnLKysrJQzPrOGWlmDrVlizpvd1aanuEVlDFQ34/xX4x8DVwNPAfxqyjoYsy52iMyIORMT2iNi+adPQVlKz7lhagj174PhxiOh93bPHQd9mUijgR8SPIuKViDgD/Am99M2gE8DmvvuXAyeLbNesM/btg1Onzl126lRvudmUCgV8SZf03f1V4OEhq90PXCHpDZLWA9cDdxXZrllnPPnkdMvNRpimLfNO4K+AN0o6IelDwH+U9F1JDwHvBP5Vtu6lkg4BRMRp4GbgG8AjwJ9FxPdK/j3M2mnLlumWm40w8dQKEXHDkMW356x7EtjVd/8QcF7LppmNsX9/L2ffn9ZZWOgtN5uSz7Q1S9nu3XDgACwugtT7euBAb7nZlJKePM3M6AV3B3grgY/wzargXnpLgI/wzeZttZd+NQ+/2ksPPnK3SvkI32ze3EtviXDAN5uH/hTO8ZzrUaTQS+9UU6c4pWNWtsEUTp66e+mdauocH+GblW1YCmdQCr30TjV1jgO+WdlGpWpS6qX3tA2d45SOWdm2bBmet19chCeeqHw4ufLGWXeqyebGR/hmZdu/v5ey6ZdCCmdQU8ZppXHANytbU6ZDqGKc7gJKiiJyr0VSu+3bt8fy8nLdwzCzWQzrVlpYSHPn1yKSjkTE9mGP+QjfzObDXUDJccA3s/lwF1ByHPDNbD588ZbkOOCbtV1dhVN3ASXHAd+szVYLp8ePQ8TZ6ROqCPrz7AJy989MHPDN2qzuwunu3b2TzT7/+d7997+/eICucyfWcBO3ZUq6A7gWeCYirsqWfRz4ZeBl4DHgAxHxwpCffQJ4CXgFOJ3XMjTIbZlmBa1Z0wuKgyQ4c6aaMZTdnrl1azPOZK5JWW2ZnwF2Diw7DFwVEW8Gvg98bMTPvzMirp402JvZgFnSGCkUTsv+lOHun5lNHPAj4l7g+YFl90TE6ezufcDlJY7NzFbddFMvHTJtGiOFwmnZATqFnVhDlZnD/yBwd85jAdwj6YikPaOeRNIeScuSlldWVkocnllDLS3Bbbedn5qZ5Cg5hWkeyg7QKezEmioiJr4BW4GHhyzfB3yZrCYw5PFLs68/DRwF3jHJ9rZt2xbWQAcPRiwuRki9rwcP1j2iZltcjOiF+/NvUt2jG+/gwYiFhXPHvbBQ7H3h91guYDlyYmrhI3xJN9Ir5u7ONjZsp3Iy+/pMtmPYUXS7lih3UJRvVOpjw4b02xPn8SljtfvnzJneV8/NM5GpJk+TtBX4apzt0tkJfAL4hYgYmn+RdAGwJiJeyr4/DNwaEV8ftz136TSQOyjKl/eaAqxfDy+/fPa+JyfrvFK6dCTdCfwV8EZJJyR9CPg0cBFwWNKDkm7L1r1U0qHsR18HfEvSUeA7wNcmCfbWUO6gKN+wnLUEF154brAHT05mI03TpXNDRFwSEesi4vKIuD0ifiYiNkev3fLqiNibrXsyInZl3z8eEW/Jbm+KCFdW2qxpHRRNOGNzWErk85+Hv//74et752o5fKatlatJHRRNqjcMy1k3bedqtXPAt3JV2Qa4tAQXX9zbjtT7fppgXfe0A0WltHMd/KR0003pf3Lqorz2nRRubsu0XAcPRqxbd36b4vr1k7foSc1tdVyVQnvisLbLwVvRNkybGCPaMn2JQ2umUZ0rk3YEuaOoHKP+Fv38ulbClzi05pi0iDqqMDlp0TKllEiTTfp6u5hcOwd8S8c0RdRRhclJi5YpTDvQBpO+3i4m184B39IxTRF1/35Yt+785evXT3eE7jM2ixv2SWmQPzklwQHf0jHNSVu7d8Of/ils3Hh22caNcMcdZ4N2yj32RceW0u827JPShz/sT04pyqvmpnBzl07H5E0Strg4/XPNY8KushQd27CfX+048kRinYe7dKwRyrwyUsodOEXHNq4rxvPpdJq7dKwZyiyipjynT9GxjVuvSSePWaUc8C0tZRVR5zHtQFl586Jjm2S9UTuFlPL/VikHfEtD2UGo7B77MufdKTq2Sbpi8nYKTZo/yMqXl9xP4eaibUfMq8Ba5rQDZRaUyxjb6s/3F2wnee3K/j0sObhoa0lLucC6as2a868pC71aw5kzo392aamXU3/yyd6R9/795RZUp3n+Ir+HNYKLtpa2Ogqs06aQNmwYvnxcPr2KFMo0dY+2TKnsOsRMHPCtflUHoWmD8NISvPji+csnOas3tSmY2zB/kOsQs8vL9aRwcw6/I6o+SWraPHbe+hs3jt9WilMwpzClchGuQ4zEiBz+NNe0vUPSM5Ie7lu2QdJhST/Ivr4252d3SnpU0jFJt5Swn7I2qXoSs2lTSHnLn39+/LZSTKE0ff6glM+xSNw0KZ3PADsHlt0C/EVEXAH8RXb/HJLWAn8EvBe4ErhB0pUzjdbaq8ogNG0QLhK025BCSU2KO9GGmOYi5vcCg4c01wGfzb7/LPArQ350B3Asehczfxn4QvZzZvWYNggXCdpVfXrpUhHTO9HZ5eV6ht2ArcDDffdfGHj8b4f8zK8D/63v/vuBT0+yPefwbW6mzWNXkfeedRspTxQ3L02vQ8wRZfXhS9oKfDUirsruvxARr+l7/G8j4rUDP/MbwHsi4l9k998P7IiI38rZxh5gD8CWLVu2HZ/k0mlmTVdk4rgmnMdglZlnH/6PJF2SbeQS4Jkh65wANvfdvxw4mfeEEXEgIrZHxPZNmzYVHJ5ZQxRp38wrVh4/Xl1qZ5KUUpfSTokqGvDvAm7Mvr8R+MqQde4HrpD0Bknrgeuzn7MU+Z+yHqOC9ri/xahiZRX96ZP0xbt3Pg15uZ7BG3An8DTwE3pH7R8CNtLrzvlB9nVDtu6lwKG+n90FfB94DNg36Tadw69YF3PBqcjrLZ/kbzHs71Zlf/okffHuna8MnkvHJuJccH2G5fCHGfxbrM6jM6rWNe95ciaZn0eqZ2wd5Ll0bDI+oaU+g+2befr/Fv1pklHm3Z8+ri9+aSn/d3LvfKUc8O2svH++DRuanddvSl2i/+SzxcXh6/T/jYYVegdV0Z8+ri9+3778TwDuna9WXq4nhZtz+BUblgtevz5i3brm5vVTqUvM0vc/btx58/SsztVTZX/6qN9v1DitdIzI4dce1EfdHPBrMPiPu3Fjs4ttKRQLZ93pjNtJpPC7TaIp42yJUQHfRVsbrekXzEhh/PMqhhc5WatKTRlnS7hoa7Nr+kRVKYx/XsXwqmcZnVVTxtkBDvg2WtMnqkph/JPudGYpLjdlquOmjLPlHPBttKYfnaUw/kl2Oj4T1SrgHL5ZFcZdaNwnvVlJnMM3q9u4lIZPejPmf8qIA75ZCkoqLjflHDM7XxVZPQd8sxSUUFx2GaDZisyQPSkHfLMUlFBcriJg2PxUkdVzwDcbpcocScHWRZcBmq2KU0Yc8M3yNCxHksI5Zja7Kk4ZccA3y9OwHEkK55jZ7Ko4ZcR9+GZ5UpiHZ0rj2v2t/dyHbzaLBuZI2j6DgdtOi3HAN8vjHElPIlG2YSWVJBUO+JLeKOnBvtuLkn57YJ1rJP1d3zq/V3S7ZnOXwjw8dUsoyjaspJKkUnP4ktYCPwR+LiKO9y2/BvidiLh2mudzDt+sZgnN8dPAkkotqszhvwt4rD/Ym1mDJdTc38CSSnLKDvjXA3fmPPZ2SUcl3S3pTXlPIGmPpGVJyysrKyUPz8ymklCUdUmluNICvqT1wPuA/zHk4QeAxYh4C/CHwP/Me56IOBAR2yNi+6ZNm8oanqUmkUKgjZFQlHVJpbgyj/DfCzwQET8afCAiXoyIH2ffHwLWSbq4xG1bkyRUCGyiSveViUXZtredzluZAf8GctI5kl4vSdn3O7LtPlfitq1J3G4xs3nsK8fuQBxlW6OUgC9pAfhF4Et9y/ZK2pvd/XXgYUlHgU8B10fKp/gmqFUZkIQKgXWa5W9a9r7SH7a6xVMrNMDqP2X/P/rCQoPzlwm1+tVl1r9p2a2J/lO0j6dWaLjWZUASKgTWZda/adlNM/6w1S0O+A3Qun/KxAqBdZj1b1r2vjKhrkurgAN+A7Tyn7KqQmCixY9Z/6Zl7yv9YWs2ib6txouIZG/btm0Lizh4MGJhIaKXve3dFhZ6y22EhF+4lIZ28GDE4mKE1PuawMuTtJT+dsMAy5ETU2sP6qNuDvhn+Z9yBouL5/5Xrt4WF+seWUT4b9pUib+tRgZ8d+lYe3m2LZuD1N9W7tKxbkqk+NHYfG+fNvwOZUnkbTUTB3wrXyrRIYGKZBtObGrD71CmBN5Ws8vL9aRwcw6/gVKraNWcKJ8l35tabj/1nHUdUvsb9cM5fKuMT908x7T53hTPqk49Z23ncg7fqtPAs8TmmYGaNt+b4lnVGzYMX96EnLWdywHfytWwita889PT5ntT218uLcGLL56/fP36huSs7RwO+FauhlW05n1EPe2ZsantL/ftg5/85PzlF13UqZkwWsMB38rVsHlyqjiinmYWidT2l3mvw/PPVzsOK4cDvpWvQRfMSO2IOrX9ZUqvTyrdvk3mgG+dltoRNaS1v0zl9fG5AOVwwLdOS+2IOjWpvD4pdi81kfvwzSx5Phdgcu7DN7PkTJOTL1JLcO7/rLIuYv6EpO9KelDSeYfk6vmUpGOSHpL01jK2a2b1myWgTpuTn7WW4Nz/gLw5F6a5AU8AF494fBdwNyDgbcC3J3lez6VjlrZZp06qao6hLs4DxIi5dKpK6VwHfC4bz33AayRdUtG2zQpzWmC4WYups5z/MEv3UmpnLtetrIAfwD2SjkjaM+Txy4Cn+u6fyJadR9IeScuSlldWVkoanqWgqUHTaYF8swbUqvr7UzqPIAVlBfx/GhFvBd4L/Kakdww8riE/M7Q9KCIORMT2iNi+adOmkoZndWty0HRLYL5ZA2pV/f2pnEeQilICfkSczL4+A3wZ2DGwyglgc9/9y4GTZWzbmqHJQdNpgXyzBtSq+vtTOY8gFYX78CVdAKyJiJey7w8Dt0bE1/vW+WfAzfSKtz8HfCoiBncK53Effns0uY/aU/yPtrTU23E/+WTvyH7//u4G1BTMuw//dcC3JB0FvgN8LSK+LmmvpL3ZOoeAx4FjwJ8AN5WwXWuQKnKp86oROC0wWkpTQdhoryr6BBHxOPCWIctv6/s+gN8sui1rrv37h1/JqaygOXilqNUaARQPQKs/76NYazpPrWCVmedH/7y0y8aN8Oyz5WzDrAk8tYIlYZ4f/fMKqM8914xOILMqOOBbK4yqBRTtBGrq+QNmgxzwrRVG1QKKtE82+fwBs0EO+NYKu3fDBRcMf2zDhtmft8nnD5gNcsC31nj1q8t/Tp90ZW3igG+tkXdh7SIX3PZcLNYmDvjWGvMIzimcdOWisZXFAd9aYx7Bue65WFw0tjI54FtrlBmc+4+q9+3r7TTqmDrARWMrkwO+tUoZJ3eldFTtovH0nALL54BvNiClo+rUi8apBdeUdtYpcsA3G5DSUXUKReM8KQbXlHbWKXLANxuQ0lF13UXjUVIMrintrFPkgG82oMhR9TxSHKnON59icE1pZ50iB3yzAbMeVaeY4pinFINryimwFDjgd0BqhbUmmOWoOi/F8ZGPpPv6F3lvpBhcU06BJSEikr1t27YtrJiDByMWFiJ6x5y928JCb7mVSzr3dc67pfL6l/HeOHgwYnGx97svLqbxe3UdsBw5MbWMi5hvBj4HvB44AxyIiE8OrHMN8BXg/2aLvhQRt457bl/xqjhfgLs6ea/1MFW+/nlXGvN7o51GXfGq8DVtgdPAv4mIByRdBByRdDgi/npgvW9GxLUlbM+mkGJhra2GXbc3T1Wv/6hr/fq90T2Fc/gR8XREPJB9/xLwCHBZ0ee1cqRYWGurYfnjjRuHr1vV6z+qddLvje4ptWgraSvws8C3hzz8dklHJd0t6U0jnmOPpGVJyysrK2UOr5NSLKy12WCx95OfrPf1H3UU7/dGB+Ul96e9ARcCR4BfG/LYTwEXZt/vAn4wyXO6aFsOF9bqVefrv7g4vHC8uDj/sfl9Vw/mWbQFkLQO+CrwjYj4xATrPwFsj4hnR63noq1ZMYM5fOgdxc+7VbGu7droom3hlI4kAbcDj+QFe0mvz9ZD0o5su88V3baZje6lr6svPcVpF4xS2jJ/Hvgm8F16bZkAvwtsAYiI2yTdDHyYXkfP/wP+dUT8n3HP7SN8s9FSPZJes6aXPBok9eobNj+jjvBLSenMiwO+2Wip9tKnOq4umGtKx2yePC3EaKn20rsDKE0O+Jasrk1GNot59NKXsZP1nDZpckrHkuW0wHhl5/BTrQnY5JzSsUZKNV2RkrKPpN1d025lzKVjNhdbtgw/wvep/+favbu8o2/vZNvNR/iWrLILfy4Aj+f5ddrNAd+SVWa6wgXgybi7pt1ctLVOcAF4cnnz51sz+MQr6zyf+Wld4S4d67xUc9OuK1iVHPCtE1LMTbuuYFVzwLdOSPHMT/e8W9WcwzeriesKNg/O4ZslKNW6wqRcf2geB3yrXNsCxay/T4p1hUm5/tBQedc+TOHma9q2z8GDEQsL515fdWGhudc7Lfr7NPW6r+OulWv1Yd7XtJ0X5/Dbp20nQLXt95mU6w/pcg6/o1JMnbRtcq62/T6Tanr9oasc8Fsq1Rxr2wLFrL9PijvjaTS5/tBpebmeaW7ATuBR4Bhwy5DHBXwqe/wh4K2TPK9z+LOrIsc6S/7ZOfz2vAZNrT+0HSNy+GUE+7XAY8A/AtYDR4ErB9bZBdydBf63Ad+e5Lkd8GcnDQ/4UjnPXyRotS1QTPv7uOBp8zQq4Bcu2kp6O/AHEfGe7P7Hsk8O/75vnT8G/jIi7szuPwpcExFPj3puF21nN+9iYleLlWVwwdPmad5F28uAp/run8iWTbsOAJL2SFqWtLyyslLC8Lpp3jnWrhYry9C2OoY1RxkBX0OWDR6/TLJOb2HEgYjYHhHbN23aVHhwXTXvuWMctGbngqfVpYyAfwLY3Hf/cuDkDOtYyXbv7qVXzpzpfS1zojAHrdmlOJGbdUMZAf9+4ApJb5C0HrgeuGtgnbuAf66etwF/Ny5/b+lavSLSqVOwdm1vmYPWdOa5MzbL86qiTxARpyXdDHyDXsfOHRHxPUl7s8dvAw7R69Q5BpwCPlB0u1aP1f7+1Wl9X3nl7JG9g5ZZ2jy1gk3F3Tnl8bVjbR5GdekUPsK3bnF3TjkGPymtngkNDvo2P55awabi7pxy+GpXVgcHfJuKu3PK4U9KVgcHfJuKWwrL4U9KVgcHfJuaWwqL8yclq4MDvlkN/EnJ6uAuHbOa7N7tAG/V8hG+tU7TLy5iNi8+wrdWcX+7Wb7WHeH76K7b5tnf7veWNV2rjvB9dGfDpn2A4v3tfm9ZG7TqCN9nL3bb0lKv42WYov3tfm9ZG7Qq4PvsxW7bty//0oFF+9v93rI2aFXA99mL3ZYXfCOKp1383rI2aFXA99mL3ZYXfBcXiz+331vWBq0K+D57sdvGBeUiXTZ+b1kb+AIo1ip5FxUZ7LKB3s7AQdvaZtQFUFp1hG/NMM9+9ryJ3dxlY1awD1/Sx4FfBl4GHgM+EBEvDFnvCeAl4BXgdN7ex9qvrn52d9mYFT/CPwxcFRFvBr4PfGzEuu+MiKsd7LutriNtd9mYFQz4EXFPRJzO7t4HXF58SNZmdR1pu8vGrNwc/geBu3MeC+AeSUck7Rn1JJL2SFqWtLyyslLi8CwFdR1pu8vGbIIuHUl/Drx+yEP7IuIr2Tr7gO3Ar8WQJ5R0aUSclPTT9NJAvxUR944bnLt02sfdMmbzNapLZ2zRNiLePebJbwSuBd41LNhnz3Ey+/qMpC8DO4CxAd/ap79rZrB10szmq2iXzk7go8AvRMSpnHUuANZExEvZ978E3Fpku9ZsvtKTWT2K5vA/DVwEHJb0oKTboJfCkXQoW+d1wLckHQW+A3wtIr5ecLtmZjalol06PxMRm7N2y6sjYm+2/GRE7Mq+fzwi3pLd3hQR7ouwqfjCI2blaNUFUKx9fOERs/J4agVLmqdEMCuPA74lzVMimJXHAd+S5ikRzMrjgG9J85QIZuVxwLekeUoEs/K4S8eS5xO1zMrhI3wzs45wwDcz6wgHfDOzjnDANzPrCAd8M7OOGHsBlDpJWgGO1z2OKVwMPFv3IGbgcVeniWMGj7tqRca9GBGbhj2QdMBvGknLTbxIu8ddnSaOGTzuqs1r3E7pmJl1hAO+mVlHOOCX60DdA5iRx12dJo4ZPO6qzWXczuGbmXWEj/DNzDrCAd/MrCMc8Esm6d9JekjSg5LukXRp3WOahKSPS/qbbOxflvSausc0jqTfkPQ9SWckJd96J2mnpEclHZN0S93jmYSkOyQ9I+nhuscyDUmbJf1vSY9k75GP1D2mcSS9WtJ3JB3NxvxvS9+Gc/jlkvRTEfFi9v2/BK6MiL01D2ssSb8E/K+IOC3pPwBExEdrHtZIkv4JcAb4Y+B3ImK55iHlkrQW+D7wi8AJ4H7ghoj461oHNoakdwA/Bj4XEVfVPZ5JSboEuCQiHpB0EXAE+JWUX29JAi6IiB9LWgd8C/hIRNxX1jZ8hF+y1WCfuQBoxB41Iu6JiNPZ3fuAy+sczyQi4pGIeLTucUxoB3AsIh6PiJeBLwDX1TymsSLiXuD5uscxrYh4OiIeyL5/CXgEuKzeUY0WPT/O7q7LbqXGDwf8OZC0X9JTwG7g9+oezww+CNxd9yBa5jLgqb77J0g8ALWFpK3AzwLfrnkoY0laK+lB4BngcESUOmYH/BlI+nNJDw+5XQcQEfsiYjOwBNxc72jPGjfubJ19wGl6Y6/dJGNuCA1Z1ohPf00m6ULgi8BvD3z6TlJEvBIRV9P7hL1DUqlpNF/icAYR8e4JV/3vwNeA35/jcCY2btySbgSuBd4ViRR3pnitU3cC2Nx3/3LgZE1j6YQsD/5FYCkivlT3eKYRES9I+ktgJ1BawdxH+CWTdEXf3fcBf1PXWKYhaSfwUeB9EXGq7vG00P3AFZLeIGk9cD1wV81jaq2sAHo78EhEfKLu8UxC0qbV7jhJ/xB4NyXHD3fplEzSF4E30useOQ7sjYgf1juq8SQdA/4B8Fy26L7Uu4sk/Srwh8Am4AXgwYh4T62DGkHSLuC/AGuBOyJif70jGk/SncA19Kbr/RHw+xFxe62DmoCknwe+CXyX3v8iwO9GxKH6RjWapDcDn6X3/lgD/FlE3FrqNhzwzcy6wSkdM7OOcMA3M+sIB3wzs45wwDcz6wgHfDOzjnDANzPrCAd8M7OO+P/djguTyiWd7wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAR/UlEQVR4nO3df4wc513H8c/n9nypLg40+K5VsH13LnIBg5qSLGkRv4Lstnb6h6lUUNyobSKQiROX8l9SKmikqlJBFFVV0pqjtdJyRy2kpq2JXMIPASmqWnyu0iRO5HA4sX11FF8IAlIjBcdf/ti93Nzc7O6sb+zdffp+SZu7mXn2me8+88wn47nbW0eEAACDb6jXBQAAqkGgA0AiCHQASASBDgCJINABIBHDvdrx2NhYTE1N9Wr3ADCQjh079mJEjBdt61mgT01NaW5urle7B4CBZPtUq23ccgGARBDoAJAIAh0AEkGgA0AiCHQASETHQLd90PY520+22G7bn7E9b/tx2zdUXyZ6bXZWmpqShoYaX2dne10R+haTpWfKXKE/KGlnm+27JG1tPvZK+tzay0I/mZ2V9u6VTp2SIhpf9+7lPEUBJktPdQz0iHhU0kttmuyW9KVo+Lak19u+rqoC0Xsf/ah0/vzKdefPN9YDKzBZeqqKe+gbJZ3JLC80161ie6/tOdtzi4uLFewaV8Lp092txw8xJktPVRHoLlhX+KkZETEdEfWIqI+PF75zFX1oYqK79fghxmTpqSoCfUHS5szyJklnK+gXfeITn5BGR1euGx1trAdWYLL0VBWBfljSB5q/7fJ2Sf8VEc9X0C/6xG23SdPT0uSkZDe+Tk831gMrMFl6yp0+U9T2lyXdLGlM0guSPiZpnSRFxAHblnS/Gr8Jc17SHRHR8a9u1ev14I9zAUB3bB+LiHrRto5/bTEi9nTYHpLuvsTaAAAV4Z2iAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BEEOgAkAgCHQASQaADQCIIdABIBIEOAIkg0AEgEQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkolSg295p+4Ttedv3Fmz/Udt/bft7to/bvqP6UgEA7XQMdNs1SQ9I2iVpm6Q9trflmt0t6amIuF7SzZI+ZXuk4loBAG2UuUK/SdJ8RJyMiFckHZK0O9cmJF1j25LWS3pJ0oVKKwUAtFUm0DdKOpNZXmiuy7pf0k9LOivpCUkfjoiL+Y5s77U9Z3tucXHxEksGABQpE+guWBe55XdJekzSj0t6q6T7bf/IqidFTEdEPSLq4+PjXZYKAGinTKAvSNqcWd6kxpV41h2SHoqGeUnPSvqpakoEAJRRJtCPStpqe0vzB523Sjqca3Na0nZJsv1GST8p6WSVhQIA2hvu1CAiLtjeL+kRSTVJByPiuO07m9sPSPq4pAdtP6HGLZp7IuLFy1g3ACCnY6BLUkQckXQkt+5A5vuzkt5ZbWkAgG7wTlEASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BEEOgAkAgCHQASQaADQCIIdABIBIEOAIkg0AEgEQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQiFKBbnun7RO2523f26LNzbYfs33c9j9XWyYAoJPhTg1s1yQ9IOkdkhYkHbV9OCKeyrR5vaTPStoZEadtv+Ey1QsAaKHMFfpNkuYj4mREvCLpkKTduTbvk/RQRJyWpIg4V22ZAIBOygT6RklnMssLzXVZb5Z0re1/sn3M9geKOrK91/ac7bnFxcVLqxgAUKhMoLtgXeSWhyXdKOndkt4l6Q9sv3nVkyKmI6IeEfXx8fGuiwUAtNbxHroaV+SbM8ubJJ0taPNiRPxA0g9sPyrpeknPVFIlAKCjMlfoRyVttb3F9oikWyUdzrX5uqRftj1se1TS2yQ9XW2pAIB2Ol6hR8QF2/slPSKpJulgRBy3fWdz+4GIeNr230h6XNJFSZ+PiCcvZ+EAgJUckb8dfmXU6/WYm5vryb4BYFDZPhYR9aJtvFMUABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BEEOgAkAgCHQASQaADQCIIdABIBIEOAIkg0AEgEQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BElAp02zttn7A9b/veNu1+3vartt9bXYkAgDI6BrrtmqQHJO2StE3SHtvbWrT7I0mPVF0kAKCzMlfoN0maj4iTEfGKpEOSdhe0+5Ckr0g6V2F9AICSygT6RklnMssLzXWvsb1R0nskHWjXke29tudszy0uLnZbKwCgjTKB7oJ1kVv+tKR7IuLVdh1FxHRE1COiPj4+XrJEAEAZwyXaLEjanFneJOlsrk1d0iHbkjQm6RbbFyLia1UUCQDorEygH5W01fYWSd+XdKuk92UbRMSWpe9tPyjpYcIcAK6sjoEeERds71fjt1dqkg5GxHHbdza3t71vDgC4MspcoSsijkg6kltXGOQRcfvaywIAdIt3igJAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BEEOgAkAgCHQASQaADQCIIdABIBIEOAIkg0AEgEQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BElAp02zttn7A9b/vegu232X68+fiW7eurLxUA0E7HQLddk/SApF2StknaY3tbrtmzkn41It4i6eOSpqsuFADQXpkr9JskzUfEyYh4RdIhSbuzDSLiWxHxn83Fb0vaVG2ZAIBOygT6RklnMssLzXWt/JakbxRtsL3X9pztucXFxfJVAgA6KhPoLlgXhQ3tX1Mj0O8p2h4R0xFRj4j6+Ph4+SoBAB0Nl2izIGlzZnmTpLP5RrbfIunzknZFxH9UUx4AoKwyV+hHJW21vcX2iKRbJR3ONrA9IekhSe+PiGeqLxMA0EnHK/SIuGB7v6RHJNUkHYyI47bvbG4/IOkPJW2Q9FnbknQhIuqXr2wAQJ4jCm+HX3b1ej3m5uZ6sm8AGFS2j7W6YOadogCQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJIJAB4BEEOgAkAgCHQASQaADQCIIdABIBIEOAIkg0AEgEQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASASBDgCJINABIBEEOgAkgkAHgEQQ6ACQCAIdABJBoANAIgh0AEgEgQ4AiSDQASARBDoAJKJUoNveafuE7Xnb9xZst+3PNLc/bvuG6kuVZmelqSlpaKjxdXa2u/Z33bW8PDbWeLTads01kt14jI01+iqz/25qXGprN9pn97djhzQ83FgeHm7Ul3/u+vXLzxkakl73uuXlWm35+2ybVv1l3XVX8b6z9Wa328WvNdtPrSZddVX7Mc0fk9nZ1bXs2NFos9TPNdesXG73+mZnV7YdG2t/3PO1ZjuaHftdjXlRdsiO1W1aTYjsDjODVjhvilZmB6To0WJgZ6+6Q1N+TkO+qKnaGc3uOLj6YGYPQNHBaDe5Z2eLBy87+fIHb2li7NixupbsPvPPyZ8g2ZNnaZJkTw67sdzqZCx70rZqV/IYd7WvSxURbR+SapL+XdKbJI1I+p6kbbk2t0j6hiRLeruk73Tq98Ybb4xuzMxEjI5GSMuP0dHG+rLt1/Ko1SLWrWu//25qvJT69u1bfu7Q0Npf01J/Wfv2Fbfdvr1zvdnX2qqf7GN4ePWY5se8qtc3MxMxMnLpfa1b13xtMzMxs+72WKf/XdVmZCR3rMsc5NHRmNn3zdXzZuT/Ymbd7asHrEyxuYGd0Z4Y1csr+9fLMaM93Q1AfgCXDvjMTPnaunmMjKx9EuQnVP5kLHvStmq3b1/5E6PbEGtB0lxEi7xuteG1BtIvSHoks/wRSR/JtfkzSXsyyyckXdeu324DfXKyeKwmJ7trX/Uju/9uaryU+mq1al/bUn9Zaz1/ll5rledhFa+vijGbnGx0NKlnS82HsjudrJ0p7kvPVjIQreqtpP/JySt3slXxyJ+MZU/aVu3KTvR249QqxFpoF+hubG/N9nsl7YyI324uv1/S2yJif6bNw5I+GRH/0lz+B0n3RMRcrq+9kvZK0sTExI2nTp0q/S+JoaHGq19dn3TxYvn2Vcvuv5saL7W+iGpfW74fe239Lb3WtfZTlaXXV8WY2dJFDWkoLiha3K1ccaxL7nRIrxb2Z13URdXWUvLl73/pQF+Jk60K+ZOx7Em71gnUbpxahVjLrnwsIupF28rcQy86NfNVlWmjiJiOiHpE1MfHx0vsetnERDXrq5bdTze1XEp9tdqlP7ddf53WdWOptrX2U4VsDVWM2cRE4z8TOt2+TZc7naidLV7fZj/daNVPJf1PTFy5k60K+VrLnrSt2pWd6O3Gqcrxa3XpvvRQn9xy4R4699DX8vq4h849dO6hR0jSsKSTkrZo+YeiP5Nr826t/KHov3bqt9tAXxrTyckIu/G10zjk2+/bt7y8YUPj0Wrb+vXLY75hw/Lx6LT/bmpcais12mf3t3378lyu1VaH78xMxNVXLz/HjrjqquXlosBf2kdRf1n79hXvO1tvdrtU/Fqz/QwNrcyDojHNH5OZmdW1bN/eaLPUz/r1K5fbvb6ZmZVtN2xof9zztWY7mtnwodigcyFdDOni6jatJkR2h5lBK5w3RSuzA1L0aDGwMyO3x6SeDevVmBw6HTPbv7D6YGYPQNHBaDe5Z2aKBy87+fIHb2libN++upbsPvPPyZ8g2ZNnaZJkTw6psdzu6q/MSduqXclj3NW+2mgX6B3voUuS7VskfVqN33g5GBGfsH1n8wr/gG1Lul/STknnJd0RufvnefV6Pebm2jYBAOS0u4c+XKaDiDgi6Uhu3YHM9yHp7rUUCQBYG94pCgCJINABIBEEOgAkgkAHgESU+i2Xy7Jje1FSu7eKjkl68QqVczlQf+8Mcu0S9fdav9c/GRGF78zsWaB3Ynuu1a/mDALq751Brl2i/l4b5Pq55QIAiSDQASAR/Rzo070uYI2ov3cGuXaJ+nttYOvv23voAIDu9PMVOgCgCwQ6ACSirwPd9n22v2/7sebjll7X1EmnD9Tud7afs/1Ec7z7/s9h2j5o+5ztJzPrfsz239n+t+bXa3tZYzst6h+IeW97s+1/tP207eO2P9xcPxDj36b+gRj/In19D932fZJejog/6XUtZdiuSXpG0jskLUg6qsYHfzzV08K6YPs5SfWI6Oc3VrzG9q9IelnSlyLiZ5vr/ljSSxHxyeb/VK+NiHt6WWcrLeq/TwMw721fp8YH2XzX9jWSjkn6dUm3awDGv039v6kBGP8ifX2FPoBukjQfEScj4hVJhyTt7nFNSYuIRyW9lFu9W9IXm99/UY2TtC+1qH8gRMTzEfHd5vf/I+lpSRs1IOPfpv6BNQiBvt/2481/mvblP90yNko6k1le0OBNkJD0t7aPNT/UexC9MSKelxonraQ39LieSzFI8162pyT9nKTvaADHP1e/NGDjv6TngW77720/WfDYLelzkn5C0lslPS/pU72stYRSH5bd534xIm6QtEvS3c1bAriyBmre214v6SuSfi8i/rvX9XSroP6BGv+sUp9YdDlFxI4y7Wz/uaSHL3M5a7UgaXNmeZOk4o9071MRcbb59Zztr6pxG+nR3lbVtRdsXxcRzzfvk57rdUHdiIgXlr7v93lve50aYTgbEQ81Vw/M+BfVP0jjn9fzK/R2mpNhyXskPdmqbZ84Kmmr7S22RyTdKulwj2sqzfbVzR8OyfbVkt6p/h/zIoclfbD5/Qclfb2HtXRtUOZ987OEvyDp6Yj408ymgRj/VvUPyvgX6fffcvkLNf7ZE5Kek/Q7S/fm+lXRB2r3tqLybL9J0lebi8OS/rLf67f9ZUk3q/EnT1+Q9DFJX5P0V5ImJJ2W9BsR0Zc/eGxR/80agHlv+5ckfVPSE5IuNlf/vhr3oft+/NvUv0cDMP5F+jrQAQDl9fUtFwBAeQQ6ACSCQAeARBDoAJAIAh0AEkGgA0AiCHQASMT/A2dPvrNrGm1mAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "\n",
    "def draw():\n",
    "    plt.scatter(x0[:, 0], x0[:, 1], c='red')\n",
    "    plt.scatter(x1[:, 0], x1[:, 1], c='blue')\n",
    "\n",
    "    plt.show()\n",
    "\n",
    "    plt.scatter(p0, np.zeros(len(p0)), c='red')\n",
    "    plt.scatter(p1, np.zeros(len(p1)), c='blue')\n",
    "\n",
    "    plt.scatter(p0.mean(), 1, c='red')\n",
    "    plt.scatter(p1.mean(), 1, c='blue')\n",
    "\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "draw()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-1"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#预测函数\n",
    "def predict(p):\n",
    "    d0 = np.power(p - p0.mean(), 2)\n",
    "    d1 = np.power(p - p1.mean(), 2)\n",
    "    pred = 1 if d0 > d1 else -1\n",
    "    return pred\n",
    "\n",
    "\n",
    "predict(p0[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.97"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#测试\n",
    "correct = 0\n",
    "for i in p0:\n",
    "    if predict(i) == -1:\n",
    "        correct += 1\n",
    "\n",
    "for i in p1:\n",
    "    if predict(i) == 1:\n",
    "        correct += 1\n",
    "\n",
    "correct / (len(p0) + len(p1))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
